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Abstract 

Results of calculations of the flow structure in binary systems with 
mass exchange driven by stellar wind are presented. 2D simulations 
have been carried out using the Roe-Oscher scheme. The fine grid we 
used allowed us to detect the details of the flow structure. In partic- 
ular, it was found out that solutions with the wind velocity V of the 
order of the orbital velocity of the system Vg are unstable: in solutions 
with V < Vs ~ the steady accretion disk takes place and if V > Vs 
- the cone shock forms. Consideration of the minor variations of the 
wind velocity V around Vg shows that they can influence greatly the 
flow structure causing transition from disk accretion to the accretion 
from the flow. During the flow rearrangement period (the accretion 
disk destruction) the accretion rate increases in dozens times for a 
short time (~ 0.1 of the orbital period of the system). This change 
can result in the drop-off of the gas from the accretor that is usually 
associated with activity in these systems. 
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1 Introduction 

The overwhelming majority of stars are binaries. In wide pairs the effect 
of one component to another is negligible but if the components are close 
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enough we can't ignore the interaction between them. One of the most 
bright manifestations of this interaction is the mass exchange between the 
components. Such systems are called interacting binaries. The study of the 
flow structure in these systems is important for the consideration of their 
evolutionary status as well as for interpretation of observations. 

Surely, the strong interaction is observed for close systems where at least 
one of components fills its critical surface (so called Roche lobe), i.e. equipo- 
tential surface passing through the inner Lagrangian point Li (see, e.g., 
In the Li point the pressure gradient is not balanced by the gravitational 
force and the component that fills its Roche lobe loses matter through the 
Li point. But in many cases even if no one of components fills its Roche lobes 
the mass transfer between them may be rather intense due to the stellar wind 
from one of the components (or in some cases from both of them). In this 
work we consider the flow structure in a binary system with mass exchange 
driven by stellar wind using symbiotic stars as the example. 

Now it is generally recognized that classic symbiotic stars are binary 
systems consisting of a cool giant and a hot compact star (white dwarf or 
subdwarf) surrounded by nebulosity. Symbiotic stars are also characterized 
by a significant photometric and spectroscopic variability and, in particular, 
by the irregular flare activity. It is also avowed that components in these 
systems don't fill their Roche lobes and there are evidences that the giant 
loses mass at the high rate due to stellar wind and the intense interaction 
between the components takes place. These hypotheses are also supported 
by observations (see, e.g., 0). 

In order to explain the flare activity in symbiotic stars different mecha- 
nisms were proposed (see, e.g., ||^). In particular, the most probable scenario 
for a classic symbiotic star is the one where the additional energy release is 
due to the fluctuations of the accretion rate near the level of the station- 
ary thermonuclear burning [Q. If the accretion rate exceeds the maximum 
value when the hydrogen can burn in a layer source of a degenerate core, 
the accreting gas stores above the burning layer and grows to the giant's size 
[p|-pT|. This process is usually associated with outbursts in classical symbi- 
otic stars. What are these changes in the accretion rate caused by? In a 
number of works (e.g., PJT^) some attempts to explain such a flare activity 
by the giant's wind variability were made but the observations carried out 
before the outburst registered only insignificant changes of the wind ||T2| , p!3 |. 
This fact imposes restrictions to all outburst scenarios caused by changes of 
giant's characteristics. 
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In this work we have carried out the 2D numerical investigation of mass 
transfer in classical symbiotic systems using TVD method that is considered 
to be the best for numerical simulation of mass transfer in binary systems |jl|] . 
The main attention was paid to consideration of the solution where stellar 
wind velocity has changed from below to above of the orbital velocity of the 
system as well as to study of the transition from the disk accretion to the 
accretion from the flow and appropriate rearrangement of the flow structure. 



2 Model 

All the calculations were carried out for a binary system with parameters of 
the classic symbiotic star - Z And. 

In order to study the gas flow structure in the equatorial plane of a 
classical symbiotic system the 2D numerical simulation has been carried out. 
The zero point of the coordinate system was placed to the center of the donor, 
X axis was directed along the line connecting centers of the components, y 
axis - along the accretor's orbital motion. The flow was described by the 
system of Euler's equations in the coordinate frame rotating with the angular 
velocity of the binary system Q: 

dp ^ dpu ^ dpv Q 
dt dx dy 



dpu dipu^ + p) dpuv 9$ 

-ir + + ^ = -Pit + ^^""p ■ 

dt dx dy dx 



dpv ^ dpuv ^ di^pv"^ + p) 9$ 2^ 

dt dx dy dy 



dpE dpuh dpvh d^ 9$ 

dt dx dy ^ dx ^ dy 



Here u = {u, v) is the velocity vector, p - pressure, p denotes density, h 
e + p/p + \uf j2 - specific total enthalpy, E = e ^ \uf jl - specific tot 
energy, e - specific intrinsic energy , $(r) - force potential. 
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In the standard problem definition when only gravitational forces from 
the point mass components and centrifugal force are taken into account the 
force potential looks as the following: 

$ r = - - ^2^ (r - . 

Ir — ril \r — 

Here Mi is the donor's mass, M2 - mass of the accretor, ri, r2 are the 
radius-vectors of the centers of components, Vc - radius-vector of the mass 
center of the system. This is so-called Roche potential. But in our case 
the additional force responsible for donor's wind acceleration should be also 
taken into account. Therefore, the form of the potential will change. 



Previous studies (e.g., [0-|T6[) have shown that the general flow structure 
in the system where components do not fill their Roche lobes is defined first 
and foremost by the stellar wind parameters. Unfortunately, wind velocity 
regime is not well-known due to the absence of an avowed mechanism of gas 
acceleration in stellar atmospheres, so we used the parametric representation 
for the force responsible for the acceleration of the wind in the following form: 

GMi r 



Fir) = a- 



where a - parameter. 

After taking into account this force directed from the donor, the modified 
force potential looks as the following: 

^, , GMi , GMi GM2 1, „2, ^2 

$ r = + a-. r - -. r - ^ 2^ {r - r^y . 

\r — ri\ \r — ri\ \r — 

Results of numerical modelling carried out for a wide range of a (see |p,7|| ) 
have shown that if a < 0.8 there is no stationary solution and gas injected to 
the space between the components falls back to the donor's surface. Finally 
there are no matter left to form the circumbinary envelope (the obtained 
densities do not exceed the background value). The situation is rather dif- 
ferent when a has larger values. In these cases stationary regime takes place 
after few orbital periods from the beginning of calculations. In this work we 
accept a = 1 as in |r^JT9[] , i.e. we assume that the accelerating force balances 
the donor's gravitational force. The accepted value of a provides the wind 
velocity change according to the /3-law (Lamers law) with /5 ~ 1 and the 
value at infinity Voo that are in agreement with observations. 
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To complete the system the ideal gas equation of state was used 

P = (7 - 1)P^ , 

where the ratio of the specific heats 7 was accepted to be ^3. The solution 
with 7 = corresponds to the case without energy losses and gives correct 
results only for a system where radiative losses are negligible, i.e. for an 
optically thick media. The observed value H/K of Civ 1550A doublet shows 
that in a quiescent state the circumbinary envelope in Z And is optically 
thin [Q. We don't consider all the envelope but only the area in the vicinity 
of the accretor in the equatorial plane of the system, namely the disk. As our 
estimates show, in the case of Z And the observed values of accretion rate is 

M„^^^ = 4.5 ■ IQ-^Me/year and the disk with T = 15000 K will be optically 
thick i.e. in the area we study the use of the adiabatic model is correct. The 
optical thickness r is defined as the product of the absorption coefficient k, 
density and the layer's width t = k, ■ p ■ I. In the case of disk accretion the 
ratio of the layer width where r = 1 to the disk width is the pacing factor. 
For the case of disk accretion 



2nR ■ H ■ pvr , 



where R 
velocity 



disk radius, H - disk half width, Vr - radial component of the 



a.sC.M 



R 



where Cg - sonic velocity, ass - Shakura-Suynyaev parameter, and 

where Vk is the Keplerian velocity, and finally this factor can be expressed 
as: 



T = l 



H 



2Tiass cl-R _ 2Txass 7^T ■ i?3/2 



Kf't^ accr ^ accr 

where TZ is the gas constant. For the typical disk temperatures T ~ 15000 K 
the absorption coefficient depends on density weakly and equals ^100 cm^/g 
g. Let us assume T = 15000 K, i? = 50/?© = A/IO, M = O.GMq (in this 
case Cs = 11 km/s, Vk = 49 km/s, H = IIRq) then we will obtain: 
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1T = 1 



H 



8-10" 



0.1 



-1 



100 cm^/g^ 



So, when Maccr = 4.5 ■ lO~^M0/yr the value of the optical width reaches the 
value r = 1 at less than 2% of the disk's width and we can conclude that the 
solution with 7 = ^3 is applicable for considered systems. 



Parameters of Z And we used in our calculations were taken from ||2T| and 
are summarized in the Table 1. The choice of namely Z And for our calcu- 
lations is rather reasonable because Z And is one of the most studied classic 
symbiotic stars. The studies of the energy distribution in the wide spectral 
range (from UV to IR) during the quiescent period (1978-1982) allowed 
to estimate the characteristics of Z And: the donor (cool M3.5III giant) loses 
mass at the rate of 2 x 10~^ Mq/jt; the hot compact component (accretor) 
has the temperature ~ 10^ K and accretes the small part of the gas of the 
donor's wind (about 2%). Gas of the circumbinary envelope has the electron 
density of ~ 2 ■ 10^° cm^^ and the temperature ~ 1.5 -r- 8 ■ 10^ K. Systematic 
observations of Z And allowed to estimate its characteristics during the out- 
burst and to find the features that can be easily explained in terms of the 
formation and subsequent drop-off of the optically thick shell by the accretor 
at the rate of about 250-300 km/s. 

Table 1. Parameters of Z And 



Parameter 


Value 


Mass of the donor Mi 


2Mq 


Radius of the donor Ri 


77Rq 


Mass of the accretor M2 


O.6M0 


Radius of the accretor R2 


0.07Rq 


Orbital period of the system Port 


760 days 


Separation A 


483/20 


Radius of the donor's Roche lobe 


232i?0 


Orbital velocity of the donor 


7 km/s 


Orbital velocity of the accretor 


25 km/s 



To solve the system of the 2D gasdynamic equations the TVD-type Roe 
scheme |22| with the restrictions of fluxes in the Oscher form was used. 
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The scheme is quasi-monotonic and has the third order of spatial accuracy 
and the first order of time accuracy. The uniform grids consisting of 301 x 301 
nodes and 601 x 601 nodes were used. The considered domain has the form 
of a square [—A . . . 2A] x [— ^/ 2A . . . ^/ 2A] with excluded circles of radii equal 
to the ones of components and centers in the component's centers. The 
free outflow boundary conditions {u = 0, p = 0) were accepted on the 
outer border as well as on the accretor surface. On the surface of the donor 
star the boundary conditions of the mass inflow were adopted. It should 
be mentioned that the boundary value of the density on the donor's surface 
doesn't influence the solution because the set of equations used allows scaling 
of the p (with the simultaneous scaling of p). Boundary conditions on the 
donor are following: density p(-Ri) = 1, temperature T(i?i) = 3200 K, the 
gas is injected to the system along the vector normal to the surface and its 
velocity is equal to V. The simulations were carried out in the velocity V 
range from 25 to 75 km/s. The solving of the set of equations was carried out 
by means of method of establishment from initial state with parameters po = 
10~^p{Ri); Po = 10~V(-Ri) c^(-Ri)/7! 1*0 = up to steady state configuration. 



3 Results 

The results have shown that when the wind velocity is less than V = 35 km/s 
the bow shock and the accretion disk behind it form in the system and the 
stationary gas flow regime takes place (see the top panel of Fig. |I]). We 
consider the stationary regime to be reached when the amount of gas injected 
to the system becomes equal to the one that leaves it (part of the gas is being 
accreted and part leaves the considered domain through the outer border). 
After that the total amount of the gas in the system doesn't change and the 
accretion rate remains constant. It happens after approximately four orbital 
periods after beginning of calculations. Results of calculations show that the 
greater the wind velocity the less the distance between the bow shock and 
the accretor and when the value of velocity is big enough the disk disappears 
and the shock becomes the attached one so the disk does not form (see the 
bottom panel of Fig. |I]). The accretion rate in these cases is slightly higher 
21- 

The observed values of the wind velocity in symbiotic stars are not far 
from the value ^ = 35 km/s that divides two principally different cases 
with different accretion regimes. So if at some moment the wind's velocity 
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Figure 1: Top panel: Density field and velocity vectors for the stationary case 
V = 25 km/s. The isolines of the density logarithm in the range of 10~^ to 1 (in 
the units of density at donor's surface) are shown. The distances are given in solar 
radii. Dashed lines show equipotential contours of the standard Roche potential. 
Donor's surface is shown by the half circle with the center in (0,0) point and the 
radius The star symbol marks accretor. Bottom panel: The same as in the 
upper panel for the case V = 75 km/s. Density contours in the bottom panel 
correspond to density values in the range of 0.04 to 0.55 (in the units of density 
at the donor's surface). 
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increases and exceeds 35 km/ s the the disk will be destroyed and the accretion 
regime will change. It should be mentioned that the value V = 35 km/s 
dividing different accretion regimes is also close to the orbital velocity (32 



km/s). This fact is in accordance with ||16|| . 

In order to consider the process of transition between these two regimes 
after the wind velocity increase we took the stationary solution for the case 

V = 25 km/s, increased the velocity up to Vi and continued calculations. 
Different cases with Vi values in the interval 35 km/s < Vi < 75 km/s were 
considered. 

In Figure |^ the morphology of the gas flow - density field and velocity 
vectors (shown by arrows) - for the case of the wind velocity change from 

V = 25 km/s to Vi = 40 km/s is presented. Donor's surface is shown by 
the circle with the center in (0,0) point and the radius Ri. Dashed lines 
mark the standard Roche equipotential contours. The star symbol marks 
accretor. Distances are given in the units of the solar radius. Two moments 
are presented: the first one is for 0.015Porb (top panel) and the second one is 
for 0.05Por6 (bottom panel) after the wind velocity change. We can see that 
after the wind velocity change two shocks (seen as condensations of density 
isolines) propagating from the donor's surface are formed. Being very close 
to each other at the first moments they begin to diverge gradually with time. 
Let us try to explain this phenomenon. Changing of the velocity regime on 
the surface of mass-losing star results in the collision of two codirectional 
flows: the 'new' one (moving with the radial velocity 40 km/s) and the 'old' 
one (25 km/s). This collision can be described in terms of the Riemann 
problem. Introducing a new coordinate system moving with the velocity of 
32.5 km/s in the radial direction it can be reduced to the frontal collision of 
two flows with velocities ±7.5 km/s. The solution of this Riemann problem 



is two diverging shocks moving with velocities (see, e.g., P5| ) 




Substituting 7 = ^3, U = 7.5 km/s, Cg = 5 km/s we obtain D12 = 
±4.5 km/s. Coming back to the original coordinate system we obtain two 
successive shocks with velocities Di = 37 km/s, D2 = 28 km/s. Note that 
these values correspond to the first moments of time after the switching of 
the velocity regime and later will be changed under the action of forces and 
variation of c,. 
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Figure 2: Density field and velocity vectors for the case after wind velocity change 
from V = 25 km/s to Vi = 40 km/s. Isolines of the density logarithm in the range 
of 10~^ to 1 (in the units of density at donor's surface) are shown. Distances are 
given in solar radii. Dashed lines show the equipotential contours of the stan- 
dard Roche potential. The star symbol marks accretor. Top panel: the moment 
0.034Por6 after the wind velocity change. Bottom panel: the moment 0.073Por6 
after the wind velocity change. 
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0.0015 




Figure 3: The pressure profile along the x axis in the area near the double shock 
wave for the moment 0.073Porb after the wind velocity change (shown in bottom 
panel of Fig. |2|). Axes in the figure: x-axis - distance - in solar radii, y-axis - 
pressure ~ in relative units. 



Let us also estimate the strength of these shocks via the relative pressure 
Ap 

jump — , where Ap = pi — po, and po, pi are the values of pressure before 
Po 



and after the shocks. Using formulae |25 




27 fiD^-Vf 



7 + 1 \ C; 



27 f{D,-V,r_^' 



2 



7 + 1 V '^i 

and substituting 7 = Y3, = 5 km/s, V = 25 km/s, Vi = 40 km/s, 

/ Ap\ ( Ap\ 

Di = 37 km/s, D2 = 28 km/s we obtain — = — = 6 (the equal 



Po Ji \Po 



2 



strength of two shocks is obvious due to the symmetry of the problem in an 
appropriate coordinate system). 

Figure |^ illustrates the pressure change in the area near the bow shocks. 
Here we can see that in the area between the two shocks the jump in pressure 
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is in accordance with the above estimates. 

When the first of shocks reaches the bow shock located in front of the 
accretor (after 0.092Porb in the case of Vi = 40 km/s) the process of the flow 
rearrangement near the accretor starts. The coming shocks begin to interact 
with the matter of the disk and crush it forcing the most part of the matter 
composing the disk to fall on the accretor's surface. 

The set of panels in Figure ^ presents the morphology of the gas flow in 
the vicinity of the accretor for different moments at early stages of the flow 
rearrangement. Here the area [0.8A . . . 1.2A] x [-0.2A . . . 0.2A] (or 80 x 80 
grid cells) is presented. After the time of about 0.09Porb (~70 days) after the 
wind velocity change the first of the two shocks formed reaches the vicinity 
of the accretor, namely the bow shock in front of it (Fig. ^ - left top panel). 
The first wave makes the bow shock to come closer to the accretor (see Fig. |^ 

- right top panel) but the disk still exists. Then ( ~ 14 days later) the second 
wave comes (see Fig. ^ - left bottom panel) and destroys the disk (see Fig. |^ 

- right bottom panel). 

The accretion rate during the flow rearrangement undergoes significant 
changes. The accretion rate behaviour during the process of the flow re- 
arrangement is presented in Figure ^ for different values of the final wind 
velocity. 

We can see that when the final wind velocity Vi equals 30 km/s there is 
no significant increase of the accretion rate while some small variations exist. 
These variations are probably due to disturbances of the accretion disk by the 
shock waves formed after the velocity increasing that are not strong enough 
to crush it. The situation is principally different for other cases - here the 
jump in the accretion rate takes place. For example, when the final wind 
velocity is equal to 40 km/s the accretion rate jumps in more than 14 times 
and in case of 75 km/s - in more than 16 times. This abrupt increase of the 
accretion rate is caused by the destruction of the accretion disk. 

It should be mentioned that the study of the transition period is limited 
in the framework of this model. After the accretion rate jump the accepted 
boundary conditions on the accretor change and the model used doesn't 
describe the real physical situation any more. Correspondingly, the presented 
results of calculations of the flow rearrangement period are correct only at 
first stages of the process. 

According to the estimates based on observations IQ in the case of Z And 
the donor's wind velocity is 25 — 40 km/s when the system is at quiescence. 
These values of the wind velocity are close to the critical value 35 km/s di- 
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Figure 4: Density field and velocity vectors in the area near the accretor for 
solution with wind velocity change from V = 25 km/s to Vi = 40 km/s. The 
arrow in the left upper corner corresponds to the velocity of 100 km/s. All other 
designations are the same as in Fig. |l[ Four moments are presented: left top 
panel - 0.092^0^6! right top - O.llOPorb! left bottom - O.llbPorb', right bottom - 
0.137Porb after the wind velocity change. 
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0.0 0.1 time "-^ 



Figure 5: The accretion rate diange witli time for the cases of wind velocity 
change from y = 25 to Vi =30, 40 and 75 km/s. Time is given in units of the 
orbital period (Porf)=758 days), the accretion rate is in percents of the matter that 
leaves donor's surface. 

viding different accretion regimes. If we suppose the presence of a minor 
variations of the wind velocity caused, e.g., by the typical giant's activity 
the wind velocity can, at some moment, become greater than the critical 
value. As the consequence the change of the flow structure and of the accre- 
tion regime will occur followed by the abrupt jump of the accretion rate in 
the transition period. 

4 Conclusions 

The 2D numerical investigation of the gas flow structure in the interacting 
binary systems with the mass exchange via stellar wind (and namely for their 
subclass - symbiotic systems) has been carried out using the Roe-Oscher 
scheme. The modelling has been carried out for different values of the wind 
velocity. The effect of the velocity change on the system's behavior was of 
special interest. 

It is shown that after the wind velocity increase the collision of two codi- 
rectional flows results in the formation of two shock waves that begin to 
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propagate from the donor's surface. These shocks can lead to the significant 
changes in the flow structure near the accreting component and namely to 
the accretion regime change. 

These results were used as the ground for the new mechanism providing 
the change of the accretion rate and explaining the nature of outbursts in 
symbiotic stars. The point of this mechanism is in the abrupt change of the 
flow structure near the accreting component as a result of minor fluctuations 
of the wind from the mass- losing component. These changes come out in the 
accretion regime change, namely, in transition from the disk accretion to the 
accretion from the flow when the wind velocity increases. The process of the 
flow rearrangement in such cases is followed by the jump in the accretion 
rate. 

Results of modelling of the flow structure in the classical representative 
of the studied type of binary systems Z And show that in the quiescent state 
characterized by the wind velocity from 25 to 40 km/s the steady accretion 
disk forms in the system. In the given regime the accretion efficiency is 
~1% of the donor's mass loss. If the wind velocity increases and becomes 
greater than 35 km/s the disk disappears and the cone shock wave forms. In 
stationary solutions for high- velocity regimes (35 — 75 km/s) the accretion 
rate is slightly higher and equals ~2% of the donor's mass loss. But during 
the flow rearrangement when the wind with increased velocity reaches disk 
and crushes it the accretion rate increases in dozens times! 

The analysis of the observations of symbiotic stars during the active stages 
shows that the most part of the registered manifestations can be explained 
in terms of the drop-off of the optically thick shell by the accretor [Q. The 
suggested mechanism provides the accretion rate jump and corresponding 
increase of the energy release rate large enough to provide the drop-off of the 
part of the gas from the accretor. Therefore it leads to the appearance of the 
flow features that can explain the observations. 
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